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Abstract: 

Monte Carlo simulations and variational calculations using a Gaussian ansatz are applied 
to a model consisting of a flexible linear polyelectrolyte chain as well as to an intrinsically 
stiff chain with up to 1000 charged monomers. Addition of salt is treated implicitly through 
a screened Coulomb potential for the electrostatic interactions. 

For the flexible model the electrostatic persistence length shows roughly three regimes in its 
dependence on the Debye-Hiickel screening length, As long as the salt content is low 
and k _1 is longer than the end-to-end distance, the electrostatic persistence length varies 
only slowly with Decreasing the screening length, a controversial region is entered. We 
find that the electrostatic persistence length scales as ^/Cp/k, in agreement with experiment 
on flexible polyelectrolytes, where £ p is a strength parameter measuring the electrostatic 
interactions within the polyelectrolyte. For screening lengths much shorter than the bond 
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length, the n~ l dependence becomes quadratic in the variational calculation. The simula- 
tions suffer from numerical problems in this regime, but seem to give a relationship half-way 
between linear and quadratic. 

A low temperature expansion only reproduces the first regime and a high temperature 
expansion, which treats the electrostatic interactions as a perturbation to a Gaussian chain, 
gives a quadratic dependence on the Debye length. 

For a sufficiently stiff chain, the persistence length varies quadratically with in agree- 
ment with earlier theories. 
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I. INTRODUCTION 



During the last decade or so a large number of simulations of charged polymers have ap- 
peared in the literature. Many of these have been devoted to the study of conformational 
properties in bulk solution, 1-6 while a few have focused on the dramatic effects polyelec- 
trolytes can produce in charged colloidal dispersions. 7 ' 8 In the former case the polymer 
chain is quite often characterized by a few geometric parameters like the end-to-end sepa- 
ration or the radius of gyration. The behaviour of a polyelectrolyte has also been analyzed 
in terms of its persistence length. The subject of persistence length may sometimes appear 
confusing due to the number of definitions and different ways to calculate the quantity. The 
general concept of persistence length is however the same — it is a measure of stiffness, a 
characteristic length over which a polymer chain remembers its direction. Of course, even a 
neutral polymer will have a non-zero persistence length and it has become common practice 
to divide the total persistence length into an electrostatic, l Pje , and a non-electrostatic or 
"bare" part, l P; o, 

lp,t = lp,0 + lp,e- (1) 



Odijk 9 and, independently, Skolnick and Fixman 10 (OSF) derived equivalent models to cal- 
culate the influence of electrostatic interactions on the persistence length and in particular 
the effect of screening from added salt. For some time there has been a controversy in 
the literature over conflicting results for the ionic strength dependence 11-14 and it is the 
purpose of this communication to resolve the conflict. 

Experimentally the OSF result, which suggests a quadratic dependence on the Debye length, 
k _1 , for the electrostatic persistence length, has been confirmed for intrinsically very stiff 
polyelectrolytes, like DNA 15 and poly(xylylene tetrahydrothiophenium chloride) 16 , but al- 
though it has been claimed that the theory is also valid for long, flexible chains, 17 ' 18 the 
experimental findings generally points towards a linear relationship, i.e. l P)e ~ k 1 , for 
flexible polyelectrolytes. 19-22 An objection to the OSF theory is the neglect of entropy. 
Taking this into account, variational calculations also give the linear relationship. 23 ' 11 ' 12 
Simulations usually show deviations from the OSF theory, and the exponent p in l P)e ~ k~ p 
may vary from 1.5 24 down to 0.3-0.9 13 . However, p « 1, in agreement with experiment, 
have also been reported. 3,6 

It is our conclusion that an almost linear relationship does exist for a flexible polyelectrolyte 
when the screening length is shorter than the dimensions of the chain, but longer than the 
monomer-monomer separation. For an intrinsically stiff chain, however, we find a quadratic 
dependence in agreement with the OSF theory. 



3 



II. THE PERSISTENCE LENGTH 



Prior to discussing the models and our different approaches we review different definitions 
of the persistence length and briefly discuss the OSF theory. Given a linear chain with N 
monomers, the microscopic approach to the persistence length is simply to sum the average 
scalar products of bond vectors with respect to a certain bond, i.e. to calculate the average 
sum of directional cosines, 

N-l-i N-l-i 

l P,i = ~ E ( r i- r i+k)=a ( cos ^i+k)- (2) 

a k=0 k=0 

where = Xj+i — Xj is the bond vector between monomer i and whereas Xj represents 
the coordinate of monomer i and a is an average bond length. In the following a will 
represent the root-mean-square monomer-monomer separation, 



N-l 



-. N-l 

3i£<rf>, (3) 



which, of course, is the same as the bond length for a chain with rigid bonds. Usually, 
i denotes the central bond and the cosines are calculated from the average of the scalar 
products in both directions and not just for subsequent bonds as written in Eq. (g). One 
may also calculated l Pt N/2 as an average over both directions. The average over all the 
local and uni-directional persistence lengths, l Pi i, is related to the mean-square end-to-end 
separation, 

N-l N-l N-l N-2 N-l 

rI = <( XJV - Xl ) 2 > = (j2r i -j2^) = H < r * > + 2 E E <*< • r i>> (4) 

i—1 j=l i=l i=l j=i+l 



i.e. 



N-l u2 



N-l ^ 2(N-l)a 2 



In the case of an infinite chain, all bonds are equivalent and so are the persistence lengths 
of Eqs. (g) and (g). We may in fact treat the polymer as infinite and disregard the index i 
in Eq. (|2|), if the chain is long enough to make end-effects unimportant, i.e., for most bonds 
(cos will become zero before k reaches the end and most of the chain will be described 
by the same orientational correlation function, Cj~. This function is normally expected to 
to decay exponentially with a characteristic length scale, l PtX , 

C k = (co S 6 k ) = e- ka / 1 ^. (6) 



Planning ahead, we might be tempted to fit the correlation function to an exponential 
function, e.g. by linear regression of InCfc vs. k, and simply calculate l PiX from the slope 
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in this case. This is not such a good idea, however, because we are disregarding one of the 
fitting parameters. It is a little optimistic to think that the fitting would reproduce the 
exact value, Co = 1, as the prefactor of the exponential (or the intercept In Co = in the 
linear regression). In order to take both parameters into account, it is better to integrate 

l P ,oo = a / C e- ka / l ^dk = C l p , x . (7) 
Jo 

We see that when Co = 1 we get l p>x back and at the same time we have a natural way to 
incorporate a fitted Co 7^ 1. The integral is essentially the sum in Eq. (|2|) and here we have 
the connection between the microscopic and exponential definitions. There is a difference, 
however. The sum in Eq. (||) is truncated at the end of the chain regardless if Ck has 
decayed to zero or not, while the integral in Eq. ([?]) extrapolates the result to that of an 
infinite chain. The significance of this will become evident when we look at the results. 



The correlation function C/% is discrete since we are considering a polymer with discrete 
bonds. We will get a continuous function, however, if we make the bonds or segment 
lengths infinitesimally small, a — > ds, and the number of segments infinitely large, while 
keeping the contour length, L = (N — l)a = X^i* a = ft) ds, fixed. This is the worm- 
like chain model. 25,26 In the chain with discrete bonds there is automatically an intrinsic 
stiffness, even in the absence of interactions. This results, in the well-known mean-square 
end-to-end distance for a freely jointed chain or random walk with a fixed step length a, 

N-l N-l 

R 2 ee = o?C + 2 J2 ( N ~ 1 " k ) a2 °k = ( N ~ (8) 

i=l k=l 

since there are N — 1 — k scalar products of the type (r, • r i+ k) [cf. Eq. (0)] and without 
interactions one has 

C fc = { 1 f ° r k = ° (9) 
I otherwise 

For the worm-like chain we need to specify C(s) = (cos 0(s)), where s = J Q S ds is a distance 
along the contour of the chain, in order to calculate [cf. Eq. (||)] 

rL 

R 2 ee = 2 (L - s) C(s) ds. (10) 
Jo 

Returning to the assumption of an exponential decay, i.e. C(s) = exp (—s/l PtW ), we get 

R 2 ee = 2Ll PtW -2l 2 p>w (l-e- L / 1 ^). (11) 



If R ee is known, l PtW can be obtained from 11. This method is often used to extract 



an apparent persistence length from experimental data, 27,3,21 although what is usually 
measured is the radius of gyration, Rq, defined as 



1 



N 

N ?J 



1 



N N 
i=l j=i+l 



(x 



(12) 
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where x cm is the centre of mass coordinate and Xij is the distance between atoms i and j. 
The worm-like chain expression used is 

R 2 G = ^£ds'(L-s>) |^,( S '- S )e-^- = ^-^ + 2 ^-2fc(l-e-^ 



(13) 

Rl^^fL. (14) 



or its long chain limit (L 3> l PtW ) 



For a worm-like chain with an exponentially decaying correlation function, Eqs. (|TT]) and 
( |i~3|) are equivalent. In other cases the obtained persistence lengths would hardly be iden- 
tical, but the trends are expected to be the same. For convenience we will be using Eq. 



The persistence lengths discussed above are all total persistence lengths, where the "bare" 
part, the persistence length in the absence of electrostatic interactions, depends on the 
chain (model) itself. For example, for the freely jointed chain l p $ is equal to the bond 
length [cf. Eqs. (||) and (||)] and, experimentally, a value extrapolated to infinite amounts 
of salt is often used. In order to calculate the electrostatic part, the OSF theory considers 
the electrostatic contribution to the bending energy, which is quadratic in the curvature or 
bending angle, for a very stiff worm-like chain (Z Pj o 3> K ~ ); whose electrostatic interactions 
are described by the Debye-Huckel approximation, i.e. a screened Coulomb potential 

e — KT 

u sc (r) = k B Tl B , (15) 

r 

where k B is Boltzmann's constant, T is the absolute temperature and l B = e 2 / {^Tre r eQk B T) 
is the Bjerrum length, with e the elementary charge, e r the dielectric constant of the 
solution and eo the permittivity of vacuum. The Debye screening length, is given 

by k 2 = 8ttI b NaIW00, where Na is Avogadro's number and / is the ionic strength in 
units of molar concentration (SI units are implied for the other quantities). The Debye- 
Huckel approximation is expected to be reasonable for a 1:-1 salt (and a weakly charged 
polyelectrolyte) , in which case the ionic strength is just the molar salt concentration, i.e. 
I = C s . The predicted electrostatic persistence length is 

hsF = [3y- 2 - 8y- 3 + e-y(y- 1 + by' 2 + 8y- 3 )] , (16) 



with y = nL and Z the amount of elementary charges on the chain. For large y this reduces 
to ^ 

l ' OSF= 4^A*> (17) 
where A = L/Z is the length per unit charge. If we translate this to a model with a certain 
bond length a and a charge per monomer a = Z/N we get 

l ' OSF = lJ^- (18) 
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Besides the bond interactions, the models used in the simulations of the flexible polyelec- 
trolyte chains only contain electrostatic interactions, which allows us to obtain a truly 
electrostatic persistence length. In the simulations of the stiff polyelectrolyte chains we 
also include a bond angle term. In a real polyelectrolyte other interactions, like excluded 
volume and attractive interactions, may be significant. Still, what we see here is in very 
good agreement with experiment. 



III. MODELS 



The polyelectrolyte is regarded as an infinitely diluted polyacid in a salt solution. The 
monomers form a linear chain where each monomer represents a charged site that carry 
one unit of negative charge. Three different models have been considered: Model 1 and 2 
are freely jointed chains with harmonic and rigid bonds, respectively. In model 3 the rigid 
bonds have been augmented with an angular potential restricting the orientation of adjacent 
bonds. The harmonic bonds correspond to an approximately Gaussian distribution of bond 
lengths and can be viewed as an averaging over the conformational freedom of the covalent 
structure that implicitly joins the charged sites. 

The solvent is treated as a dielectric continuum with a dielectric constant e r = 78.3, corre- 
sponding to water at room temperature, and the salt ions are implicitly represented by the 
screening parameter k. Thus, charged monomers interact via a screened Coulomb poten- 
tial, Ec, which for the flexible bond model (model 1) is augmented with a Gaussian term, 
Eq , leading to a total interaction energy 

where N is the number of monomers, r^j is the distance between monomer i and j, e is the 
electronic charge, eo is the permittivity of vacuum and a the (fractional) amount of charge 
on a monomer i. The force constant, k, is given implicitly through the input parameter 
tq = (e 2 /47re r eofe) 1 / 3 , which is the equilibrium distance for a fully charged dimer. 

The angular potential in model 3 is simply, 

E A = -k B T 1 -^ ^ cos e h (20) 
r o i=l 

where 0i is the angle between bond i and i + The angular potential is characterized by 
the length l A and increasing I a means increased stiffness in the chain. 
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IV. METHODS 



The key approach employed for studying the scaling properties of the persistence length is 
the Monte Carlo method. In this section, however, we also discuss what is expected from 
non-stochastic approaches; a variational method, high- and low-T expansions. 



A. Monte Carlo simulation 



The Monte Carlo (MC) simulations were performed with the traditional Metropolis algo- 
rithm 28 in a canonical ensemble. One polyelectrolyte chain at infinite dilution was simu- 
lated. The sampling of chain conformations was made highly efficient, compared to tradi- 
tional Monte Carlo moves, by using a pivot algorithm, which allows chain lengths of more 
than 2000 monomers. The pivot algorithm was first described by Lai, 29 and its efficiency 
for self-avoiding walks have been thoroughly discussed by Madras and Sokal. 30 Despite its 
excellent properties, it was only recently that the pivot algorithm started to gain wider 
recognition as a tool to enhance simulations of polymers and polyelectrolytes. 31 ' 24 ' 32 ' 33 

In a traditional single-move algorithm only one or a few monomers are translated or rotated 
at a time. The number of interactions that has to be calculated is of the order iV and a large 
number of attempts per monomer is needed to generate independent chain conformations. 
Similarly, in the pivot algorithm each monomer i (except the first one) is translated in 
turn but together with the rest of the chain (monomers i + 1 to N). Furthermore, the 
last part of the polymer is then rotated as a rigid body around one of the coordinate axes 
with monomer i as origin. To be specific, if the bonds are flexible (model 1) they are both 
translated and rotated, while in the case of rigid bonds (models 2 and 3) only the rotation 
is performed. The number of interactions calculated in one step is of the order iV 2 but 
independent conformations are obtained after a only a few attempted moves, on the order 
of one per monomer or N in total. The net effect is a greatly reduced simulation time for 
a given degree of precision and a computational cost that grows approximately as iV 3 . 

One expects the pivot algorithm to have its maximal efficiency for unscreened and extended 
chains, but we have recently shown its efficiency even for screened chains with less elongated 
structure. 34 The pivot algorithm provides an effective sampling of global properties like 
R ee , but it turns out to be superior also for local properties like the monomer- monomer 
separation. A more detailed investigation than the present, is needed to really map out the 
applicable range for the pivot algorithm. 
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B. Variational approach 



In the variational approach the true Boltzmann distribution P oc exp(—E/kBT) is approx- 
imated by a trial distribution Py oc exp(— Ey /hsT), having a set of free parameters. It 
yields a variational free energy 

F = {E) V - k B TS v > F, (21) 

bounded from below by the true free energy F. (E)y is the average of the true energy 
in the trial Boltzmann distribution Py and Sy is the corresponding variational entropy. 
By minimizing F with respect to the parameters of Ey, thermodynamic quantities can be 
estimated. In the case of model 1 (flexible chain), we consider an arbitrary pure Gaussian 
variational distribution, 32 corresponding to 

i?vA B T = ^G ij 1 r ! T i , (22) 

where rj denotes the zth bond vector. With the correlation matrix G optimized (numeri- 
cally) so as to minimize F, the variational estimate for the correlation matrix becomes 

(r i -T j )v = aGi j . (23) 

From this a variational estimate for the persistence length, e.g. as defined in Eq. (|5|), can 
be obtained. 



C. High temperature expansion 



For model 1 at high temperatures, the energy [Eq. (|l~9|)1 will be dominated by the bond 
term, and the interaction Eq can be treated as a perturbation. A perturbative expansion 
for the thermal average of an arbitrary observable / yields 

(/) = (/) ° " k^f ifEc) ° c + 2ik^jf {fEcEc) ° c + • ' ' ' (24) 

where ( )^ refers to connected expectation values in the unperturbed Boltzmann distribu- 
tion. 



In a continuum approximation, valid when the nearest-neighbour distance is small as com- 
pared to the range of the interaction, the first few terms in the perturbative expansion of 
the squared end-to-end distance R^, e are given by 35 (assuming full charge, a = 1) 



3 + 7/i(^)+7 2 / 2 (/u) + 



(25) 
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where 



fM = - l^- 1/2 + O(^) , / 2 (m) = - (- - f ) + 0(^ 1/2 ) (26) 

V 7T V 7T 9 / 



with the dimensionless parameters \i and 7 given by 

Note that /i is proportional to the squared unperturbed end-to-end distance in units of the 
screening length; in the parameter range we are studying, fx is large. 

Keeping only the leading large-yU contribution to each term /j yields a power series in 7, 
indicating that the expansion is valid for small 7. Unfortunately, for the parameter range 
we are studying 7 is not small, and while the first order expansion is proportional to k~ 2 , 
the next order term, proportional to k -4 , is larger and has a negative value. 

Similarly, for the average squared bond length we obtain in the continuum approximation 
the perturbative result 



2 = / 2 \ = Ijj 

u — V mm I i 



3 r /2 
3 + 4-^ + . 



(28) 



The persistence length 1 P: r is then obtained by dividing R^ e by Nr mm . In the results 
section, we have included the first order perturbative result for comparison, although a 
good agreement cannot be expected since the parameter 7 is not small. 



D. Low temperature expansion 



A similar approach can be followed for low temperatures, where the Boltzmann distribution 
is dominated by the minimal energy configuration. For small fluctuations around the ground 
state, the energy can to lowest order be considered Gaussian, with the corrections treated 
as a perturbation. The resulting correlation matrix (r< • rj) 36 can be used to calculate the 
persistence length. 



V. RESULTS AND DISCUSSION 



A. Flexible polyelectrolyte 

From the form of the Hamiltonian, being just a sum of screened Coulomb interactions 
[Eq. ([l5|) l, and the fact that the statistical mechanical weight of a conformation is given 
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by exp(—E/kBT), it is easy to see that the rigid bond model scales with the bond length 
a as the unit of length and that any scaled length, e.g. l p /a, must be a function of the 
dimensionless numbers (scaled length scales) na and the strength parameter 

c?ln 

t P = —, 29 
a 

which measures the strength of the electrostatic interactions within the polyelectrolyte 
for a given N. As a visual proof, Fig. 1 shows that l p ^ e /a is constant for a given value 
of these parameters. Note that £ p is different from the so called Manning parameter, 37 
£ = Ib/A = oIb/cl, which measures the interactions between the polyelectrolyte and mobile 
ions. In a previous paper we found for a titrating polyelectrolyte that the approximation 
of putting full charges a distance A apart is inferior to having fractional charges on the 
actual titrating sites. 38 Thus, we may conclude that £ p is the more pertinent parameter in 
this context. 



The results in Fig. 2 show two distinct regimes in the salt dependence of the electrostatic 
persistence length. For low salt concentrations, i.e. when the Debye screening length is 
longer than the dimensions of the chain (k _1 3> R ee ), lp,e is almost constant. Note that 
the symbols on the y-axis represent the values in the unscreened case (k = 0). When 
k^ 1 <C R ee the relationship between l PiC and k _1 is close to linear. There is actually a third 
regime when k _1 is less than the bond length, but we will save that till the discussion of 
the variational results. 

Experimentally, Reed et al. have found a linear dependence on both k~ 1 and a for l Pje . 21 
To us, this implies 

a 

(30) 




V8000vriVA/a 



at high salt concentrations, which is confirmed in Figs. 3 and 4. Note that l pe is predicted 
to be temperature independent. Ha and Thirumalai (HT) have derived a similar expression, 

l p ~ yj (Ipfils/A 2 )/ k for Iosf 3> ^p,o- 12 As long as 1 P) q = a the two expressions are the same, 
although l p in their case appears to be the total persistence length. However, if we add 
a repulsive potential between next-nearest neighbours, which increases l p $ (but not too 
much, see below), Eq. (^) still holds, while the HT expression with the actual intrinsic 
persistence length does not. The slope for the highly screened regime (k _1 <C R e e) i n Fig. 
4 is roughly -1.1. 

Equation (|3(]) has no N dependence. Naively one may expect excluded volume behaviour in 
the high salt regime as a result of short ranged electrostatic interactions between monomers, 
which can also be far apart in the sequence as the chain loops back on itself. For a self- 
avoiding chain we have R 2 ee ~ A^ 1,2 a 2 , 39 which from Eq. (|5|) would suggest l p>e ~ A r °' 2 for 
long chains. The slopes of the linear parts of Fig. 3, roughly 2.6, 4.3 and 7.3 for N = 80, 
320 and 1000, respectively, give an iV-exponent of 0.4. 
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FIG. 1. The logarithm of the scaled electrostatic persistence length, where l p ^ e = I p r — a 
[cf. Eq. (pi)], as a function of ln(Ka) for rigid bond chains with N = 320 and fixed £ p . Filled 
squares represent a = 24 A and a = 1 and the open symbols represent a/a = 0.3536/3 
A (triangles), 0.5/6 A (diamonds) and 0.707/12 A (circles), with the values on the y-axis 
showing the results for the salt-free case. The dot-dashed line marks a slope of -1. 
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FIG. 2. The logarithm of the scaled electrostatic persistence length, where l p ^ e = I p ^r — a 
[cf. Eq. (||)], as a function of ln(«a). The results are for rigid bond chains with N = 320 
and a = 3 A (triangles) and 12 A (circles) and fully charged monomers. The values on the 
y-axis are for the salt-free case, and are connected to the rest of the points by thinner lines. 
The dot-dashed line denotes a slope of -1. 

So far we have used the "macroscopic" persistence length, Eq. (||), but the same conclusions 
are reached with other definitions as can be seen in Table I and Fig. 5. The two persistence 
lengths obtained directly form the simulation, Z P) jv/2 an d lp,R [Eqs. @ and (||)], are almost 
identical. l p ^ is a global quantity and thus have better statistics than 1 P) n/2i which has more 
local character. This becomes noticeable at high salt concentrations where the latter starts 
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FIG. 3. The scaled electrostatic persistence length, where l p>e = 1 P: r — a, as a function of 
Y^/ko. The results are for rigid bond chains with N = 80 (circles) and 320 (squares) 
with bond lengths 3 A, 6 A, 12 A and 24 A as well as JV = 1000 with a = 3 A and 12 A 
(triangles). In order to clearly show the linear regime, only points where ft -1 < R ee /W are 
displayed. 




-4-2 2 

ln(xa/VQ 

FIG. 4. The logarithm of the scaled electrostatic persistence length, where l Pif , = l Pt n — a, 
as a function of y^/fta. The lines represent the same combinations of a and a as in 3. 
The dot-dashed line marks a slope of -1. 

to suffer from insufficient averaging. The worm-like chain persistence length, Eq. (p4]) ; also 
give the same results for k^ 1 <C R e e and > a. For screening lengths less than the bond 
length the true structure of the chain becomes important and the worm-like chain model 
breaks down and even gives negative persistence lengths. To save the worm-like chain, a 
term corresponding to the last term of Eq. (||) would have been needed in Eq. (|i"T|). 

lp,x [Eq. (H)] and Z P;00 have been obtained by fitting an exponential to the orientational 
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TABLE I. End-to-end separation and electrostatic persistence lengths from different defi- 
nitions, l Pt N/2 — o. from Eq. (|2|), l Pt R — a from Eq. (pi), l PjX — a from Eq. (|6|), l Pt00 — a from 
Eq- (0), l Pt w ~ a from Eq. (|il| ) and /osf from Eq. (p~6|) , obtained for a 320-mer. The unit 
of length is 1 A. 



a 


C s /M 


K- 1 




lp,N/2 ~ a 


lp,R — a 
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96 
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12 


0.01 


30 
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84 
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87 


80 
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481 


27 


24 
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14 


18 


1 




1 


3.0 


313 


7 


7 


213 


-11 


0.8 


0.1 




10 


0.96 


234 


4 


1 


262 


-12 


-5 


0.01 




100 


0.30 


217 


1 


0.2 


420 


-12 


-6 


0.001 



correlation function. In order to avoid end-effects we have calculated C& from the central 
bond and fitted the function only in the range A;=0-120 (iV=320). Figure 6 shows the 
correlation function together with the fitted exponential. It is clear that some information is 
lost at small values of k and that this is of greater importance at higher salt concentrations. 
In particular the initial part at Co = 1, which for the integrated persistence length l poo 
corresponds to the intrinsic persistence length, l p $ = a, is lost. Table I shows that if we did 
not subtract 1 Pi q from l Pt00 , this persistence length would be in good agreement with the 
microscopic definitions, not only at intermediate salt concentrations, but also at higher salt 
concentrations. Even more information is lost when calculating l px . Note for example that 
the fitted exponentials have a prefactor (Co) that is less than one and decreasing faster than 
C\. In Fig. 5 we see that l PiX has a behaviour of its own. There is not a simple, consistent 
k^ 1 dependence and an apparent exponent p in l PjX ~ k~ p can have a range of values less 
than one, which is exactly the result of Micka and Kremer. 13 Their conclusions are thus a 
consequence of an unfortunate choice of persistence length definition. 

The OSF persistence length also has an almost constant regime, with a salt-free limit 
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FIG. 5. The logarithm of the electrostatic persistence length as a function of the logarithm 
of na/y/^, for different definitions of l pe : 1 P) n/2 ~ a horn Eq. (]2|) (plus signs), 1 Pj r — a from 
Eq. (H) (squares), l P)X — a from Eq. (g) (crosses), i Pj00 — a from Eq. (0) (diamonds), l PtW — a 
from Eq. (|ll|) (triangles). The values are based on simulation results, either directly or 
through R ee , of 320-mers with a = 3 A. 
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FIG. 6. The orientational correlation function for a 320-mer with a = 12 A and salt 
concentrations of 0, 0.1, 1, 100 mM (top to bottom) represented by small circles. The solid 
lines are exponentials fitted in the range /c=0-120. The dot-dashed line just shows the 
upper limit of the fitting and the upper limit of the y-axis is a reminder that Cq = 1. 




Iosf = IbZ 2 /72 (see Fig. 7). In this region the value of the electrostatic persistence length 
is sensitive to whether it includes end-effects, like l P: i and 1 Pi r, or is defined through a 
correlation function that is extrapolated to infinite chain lengths, like l PtQO (and in principle 
l PjX if Co is close to one, cf. Eq. (0)). 1<jsf belongs to the latter group and is actually 
in fair agreement with l Pj00 (and l P)X ) in this regime (see Table I). The worm-like chain 
expression, Eq. (11), gives a behaviour intermediate between the two extremes and since 
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FIG. 7. OSF theory [Eq. ([if])] (no symbols) compared to the simulated electrostatic 
persistence length, l p>e = 1 Pi r — a from Eq. (]5|), (squares). The results are for 320-mers 
with a = 3 A (solid lines) and a = 12 A (dashed lines). The values on the y-axis are the 
results for the salt-free case. 



the worm-like chain is often used to obtain persistence lengths from experimental data 
this is the behaviour expected to be reported experimentally. Such a report by Forster et 
al. appears to be at odds with both OSF theory and the notion of a power law for the 
k dependence, 40 but a comparison with the results here suggests that their experimental 
system just displays the constant regime and the subsequent cross-over region. This means 
that their comparison with OSF theory is rather unfair since they are using the high-salt 
version, Eq. fll8j). 

In the second regime, Iosf of course displays the well-known quadratic dependence on the 
Debye length, in disagreement with simulations. It is interesting to note that I'osf/ 0, nas 
the correct relationship between £ p and na. It is essentially the square of the dimensionless 
form of the right-hand side of Eq. (^0[). This is a little bit tricky, however. Going from 
Eq. (|17]) to Eq. ( |l8| ) the choice of the bond length a is arbitrary, because it is always 
balanced by a. This is no longer true if we are to scale I'qsf by a - in that case we need 
the proper rigid bond model to get the correct bond length dependence, which is closely 
related to the above discussion about the difference between having an effective charge or 
an effective bond length, i.e. using £ p or the Manning parameter, £. 

Figure 8 shows that the electrostatic persistence length obtained from the variational cal- 
culations is slightly higher than the simulation results for the same model. Otherwise, the 
agreement is excellent. The calculations have been made for a flexible bond chain with 
tq = 6 A and the root-mean-square bond length is around 12 A, a little longer at low salt 
concentrations and a little shorter at high screening. As the simulation results show the 
behaviour is almost identical to that of a rigid bond chain with a = 12 A. It is a general 
conclusion that models with flexible, harmonic bonds give the same results as rigid bond 
chains if the the actual, simulated root-mean-square monomer-monomer separation is used 
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FIG. 8. Comparison between simulations and variational and perturbational calculations. 
The simulation results are for 320-mers with rigid bonds of 3 A (triangles) and 12 A (circles) 
and flexible bonds with tq = 6 A (dashed line, stars). The solid line represent the variational 
result for the latter, while the dotted line and the long dashed lines are low temperature 
and high temperature expansions, respectively, for a 512-mer with ro = 6 A. The symbols 
on the y-axis are for the salt-free case and the dot-dashed lines mark a slope of -1 and -2, 
respectively. 



to represent the bond length. A pair similar to ro = 6 A/a = 12 A is ro = 9 A/a = 24 A. 

When the Debye length is on the order of the bond length («a ~ 1) the variational result 
crosses over to a third regime where l pe is proportional to k~ 2 , i.e. we get a behaviour 
similar to the OSF prediction, which applies to stiff chains. Although the chain is highly 
flexible, it is stiff on the length scale of a bond and when k _1 is of this order, the intrinsic 
persistence length, l P; o = a, dominates the total persistence length. The OSF theory 
treats the electrostatic interactions as a perturbation to an intrinsic bending energy, which 
is quadratic in the bending angle. Analysing the high temperature expansion, Eq. (|24]) , 
which is also a perturbational calculation, we would get a k~~ 2 dependence to the lowest 
order for the entire range of salt concentrations if l P)e were plotted as a function of -\/£p/ k, 
but scaling with the bond length as in Fig. 8 the expansion result curves upwards at lower 
salt concentrations. 

The simulation results appear also to cross over into a third regime at (ridiculously) high 
salt concentrations, although the exponent of k _1 seems to be closer to 1.5 and slowly 
increasing. Seidel's large exponent, also 1.5, for k~ 1 may be due to the fact that the Debye 
length was smaller than the distance between neighbouring charges. 24 

Finally, the low temperature expansion, which corresponds to a rigid rod, gives the expected 
almost constant persistence length at low salt concentrations, but diverges when nR ee ~ 1 
and the chain cannot bend as it is supposed to. 
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FIG. 9. The electrostatic persistence length for a 320-mer with rigid bonds of a = 24 A and 
an angular potential with I a = \/2000 A (filled circles) and I a = v / 20000 A (filled triangles) . 
The solid lines are least square fits, which have the slopes -1.0 and -2.0, respectively. Also 
shown is the corresponding chain without the angular potential (open squares). 

B. Stiff polyelectrolyte 

Thus, we have shown that the OSF theory is not applicable for a flexible polyelectrolyte 
where the electrostatic persistence length is of the same order as the total persistence 
length. This should not come as a surprise, since the OSF theory was derived for a chain 
with a significant intrinsic stiffness and where the electrostatic persistence length only gives 
a minor contribution to the total persistence length. 

One way to introduce a non-electrostatic stiffness beyond the bond length is via an angular 
potential [see Eq. (po|)1. We have performed such simulations for a fully charged chain with 
320 monomers and a fixed bond length a = 24 A and found that the persistence length 
still scales like re . With increasing I a, however, the dependence changes over to re~ 2 in 
agreement with the OSF result. Figure 9 shows how a chain with I a = \/2000 A and with an 
intrinsic persistence length, l P: o = 83 A decays like re -1 , just like the flexible chain without 
the angular potential, i.e. Ia = A. Upon increasing the stiffness to Ia = \/20000 A and 
a Ipfl = 759 A the re dependence shifts to -2. The re dependence at high salt concentration 
still seems to be quadratic although considerable numerical difficulties prevent an accurate 
analysis. 

VI. CONCLUSIONS 

Simulating a flexible polyelectrolyte model that only has screened Coulomb interactions, 
besides the bond interactions, we find a linear relationship between the electrostatic per- 
sistence length and the Debye length in the controversial region where the Debye length 
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is shorter than the chain dimensions (k~ 4 -C R e e) and the interactions have become more 
local. Furthermore, we also find that the electrostatic persistence length is proportional to 
the square root of the strength parameter £ p = a 2 lB/a, which measures the internal inter- 
actions, i.e we have l p , e ~ \fip/^- This is in agreement with experiment, but contradicts 
OSF theory, which predicts l Vfi o?Ib / '(2/ta) 2 . 

The conclusions remain the same if different definitions of persistence length are used, 
although l P:X [Eq. @] may give erratic results when obtained by a simple fitting to the 
orientational correlation function. 

Variational calculations agree with simulations. A small difference is seen at very high 
salt concentrations, where the variational results predict a k~ 2 dependence as in the OSF 
theory. Here simulations give an exponent p between 1 and 2 for l p ^ e ~ k~ p . In this region, 
k^ 1 <C a, and the electrostatic interactions can be treated as a small perturbation to the 
harmonic interactions in the variational ansatz, similar to the perturbation of the quadratic 
bending energy in the OSF theory. A high temperature expansion carried to first order 
gives a quadratic dependence throughout the concentration range. 

Finally, for a stiff chain where the electrostatic persistence length only gives a minor contri- 
bution to the the total stiffness, we find that theory due to Odijk and Skolnick and Fixman 
is valid. Under these conditions the electrostatic persistence length scales like k~ 2 . 
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